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Incorporation  of  (5)  ar*d  (6)  into  Eq  (3)  yields 
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Preface 


This  document  was  originally  submitted  by  Capt  Edward  L.  V/olf 
as  his  master's  thesis  to  the  faculty  of  the  School  of 
Engineering,  Air  Force  Institute  of  Technology,  Air  University 
in  partial  fulfillment  of  the  requirements  for  the  degree  of 
I'tester  of  Science  (AFIT  GME/PH/76D-8 ) . The  AIRDIF  computer 
code  was  developed  by  Capt  Raymond  A.  Shulstad  at  the  Air  Force 
Weapons  Laboratory.  A previous  report,  AFWL-TR-76-221 , docu- 
m.ented  the  AIRDIF  theory  and  presented  the  results  of  a study 
where  AIRDIF  was  used  to  evaluate  mass  integral  scaling.  This 
report  presents  a user's  guide  for  AIRDIF. 
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These  transport  calculations  have  ceen  carried  out 
either  by  Monte  Carlo  sirciulation  or  by  the  nethod  of  nass 
integral  scaling  (MIS).  Shulstad  (Ref  1)  gives  a list  of  the 
Monte  Carlo  results  which  have  been  published  in  the  open 
literature,  and  he  also  describes  the  nethod  of  nass  integral 
scaling.  As  Shulstad  points  out,  MIS  is  essentially  a "fix- 
up" of  a one  dimensional  (homogeneous)  air  calculation. 

AIRDIF,  the  code  described  here,  is  an  "in  between" 
method.  It  is  a true  two  dimensional  calculation  but  one 
which  involves  several  simplifying  assumptions  which  permit 
its  execution  in  a fraction  of  the  time  required  for  a 
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Monte  Carlo  calculation.  Nevertheless,  as  Shulstad  demon- 
strates (Ref  1 ),  the  agreement  between  AIRDIF  and  Monte  Carlo 
is  excellent  for  altitudes  up  to  25 

The  basis  of  AIRDIF  is  a special  form  of  the  multigroup 
diffusion  equation  which  is  described  in  Section  II.  The 
simplifying  assumptions  are  those  inherent  in  the  use  of  the 
diffusion  equation:  nam.ely,  that  the  region  of  interest  Is 

large  with  respect  to  the  radiation  mean  free  path,  the 
spatial  variation  is  slowly  varying,  and  the  directional 
distribution  of  the  radiation  flux  is  linear  with  the  cosine 
of  directional  angle. 

The  AIRDIF  code  is  user  oriented  in  that  it  is  simple 
to  operate  and  does  not  require  an  extensive  knov/ledge  of 
diffusion  theory  or  computer  technology.  AIRDIF  produces 
radiation  envirorjoents  in  one  dimensional  (1-D)  infinite 
homogeneous  air  as  well  as  in  2-D  variable  density  air. 

Given  the  proper  input,  this  code  is  capable  of  producing 
correction  factors,  called  K-?actcrs,  vfhich  are  defined  to 
be  the  ^mr^  2-D  variable  density  air  dose  divided  by  the 
iimr^  1-D  MIS  dose. 

This  report  contains  five  sections  and  five  appendices. 
Section  II  contains  a brief  description  of  the  theory.  Sec- 
tion III  describes  the  code  itself,  and  Section  IV  is  a 
users 's  guide  (is,  instructions  on  how  to  use  this  code). 
Section  V discusses  limitations  inherent  in  the  code. 
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II.  Theory 


This  section  is  not  intended  to  give  detailed  derivations 
of  the  theory  behind  the  calcnlations  within  the  AIRDIF  code. 
Instead,  this  section  presents  a sunsiary  of  the  theory  pre- 
sented in  reference  1 . 


Atmospheric  Diffnsior.  Eauation 

Given  a nuclear  detonation  at  a specified  altitude  and 
the  proper  input  data,  AIHDI?  computes  the  multigroup  energy 
fluences,  doses,  and  K-Factcrs  (if  opted)  throughout  a spa- 
tial domain  about  the  so^orce.  AIRDI?  is  based  on  a special 
form  of  diffusion  theory/  which  is  s’ummarized  below.  An  ana- 
lytical first  collided  source  is  used  so  that  diffusion 
theory  is  only  applied  to  the  solution  of  the  collided  flu- 
ences. Solutions  are  thus  obtained  separately  for  the 
collided  and  luncollided  fluences  over  each  energy  group 
throughout  the  region  of  interest.  The  total  group  fluences 
are  then  obtained  by  simple  addition  of  the  collided  and  un- 
collided particle  fluences  every^vhere.  Doses  are  computed 
by  summing  over  all  groups  the  product  of  the  fluence  and 
the  dose  response  function  in  each  group. 

As  mentioned  above,  the  ijncollided  particle  fluence  is 
determined  analytically . The  equation  used  is 


,g-^  S|_|_  exp(-^’’*  RHGR) 


Virr' 


(1) 


where  is  the  sea  level  density  (1.225E-3  gm/cm^, 


r.*'  — 


F®  = the  uncollided  fluence  in  group  g' 

S®  = the  number  of  grout)  g'  oarticles  emanating 
0 

from  an  isotropic  point  source  at  r=0  and 

r = slant  range  between  the  source  and  mesh 
point  in  units  of  cm 

fT  * 

Zz  - the  total  macroscooic  cross  section  at  sea 

i 

- 1 


level  in  lonits  of  cm 
RHCR  = /q  p dr  in  units  of  g/cm^ 


(2) 


and  where  p is  the  atmospheric  density  as  a function  of 
altitude.  Note  that  "RHOR"  is  often  referred  to  as  the  mass 
range  (gm/cm^)  or  optical  density. 

The  collided  fluences  are  found  numerically  by  solving 
the  diffusion  equation  which  can  be  written  as 


7-  (D^7F*) 


• &T&  - 

'R' 


— _ q s 


g' 

g'  = 1 ^ 


(3) 


where 


= 


6 " 
T)S  = 

yg  = 
= 


h - 


ss  = 


the  energy  group  number 

the  group  g diffusion  coefficient  in  cm 

the  group  g fluence  in  particles/cm- 

the  macroscopic  removal  cross  section  for 

group  g in  cm" ‘ 

the  macroscopic  scatter  cross  section  from 
group  g ' into  group  g in  cm" ' . 
the  source  of  uncollided  particles  for  grout 
g in  particles/cm^  given  by 
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sity  (Pq=  1 . 225x1 0“ ^g/cm ^ ) . The  cross  section  at  any  altitude 
Z is  determined  within  AIRDI?  by  multiplying  the  sea  level 


value  by  the  ratio  of  the  densities,  f(z),  at  altitude  z and  sea 
level,  i . e . , 

f(z)  = p(z)/p(0)  (5) 

where  p(z)  is  the  density  at  altitude  Z as  given  by  the  1962 
U.S.  Standard  Atmosphere  (Ref  2)  which  was  incorporated  into 
AIRDIF. 

The  diffusion  coefficient  in  Eq  (3)  incorporates  a trans- 
port correction  and  is  defined  by 


DS  E 


1/{3ZS  f(; 
‘^0 


)}  = 


= D^VfCz) 


(6) 


where 

Z®  is 

the  mac 

roscopic 

transport  cross  section  at  sea 

level 

and  is 

f*  HZ*  "v  1*1 6 1* 

defined 

by 

j-g  - yB  - yS"^S  (y) 

si 


where  is  the  total  macroscopic  cross  section  and  z*"^*  is 
T SI 

the  cofficient  in  the  Legendre  expansion  of  the  angular 
dependence  of  the  scattering  cross  section  for  the  in-group 
scatter  of  group  g. 


5 


•< 


I 


Incorporation  of  Sq  (5)  and  (6)  into  Eq  (3)  yields 


Df  VF 


-Zf  f(z)  F® 

% 


z|'"Sf(z)  FS' 

g’  = l 


(8) 


Inspection  of  the  spatial  leakage  terai  V- 


{D«/f(z) }7FS 


Eq  (8)  reveals  that  the  density  gradient  (9/9z (p^/p^) ) aust 
also  be  deterained  as  a fionction  of  altitude.  A aodel  was 
developed  under  the  assumption  of  local  exponential  density 
variation  and  incorporated  into  AIRDIF. 

Diffusion  theory  as  employed  in  AIRDIF,  separates  col- 
lided and  uncollided  multigroup  fluences.  The  uncollided 
fluence  is  determined  by  analytically  evaluating  Eq  (1).  The 
collided  fluences  are  determined  numerically  by  solving  Eq 
(8).  The  total  group  fluences  are  then  found  by  a simple 
group  by  group  addition  of  the  collided  and  uncollided  flu- 
ences. Doses  are  determined  by 


Dose  = I FS(r)  Df 


(9) 


where 

F®(r)  = the  total  fluence  for  group  g at  location  r 

Df  = the  dose  response  function  to  a given  material 
H 

for  energy  group  g. 

Finally,  K-Factors  (if  opted)  are  computed  as  the  ratio  of 
the  2-D  ^-vr^  variable  density  air  dose  to  the  1-D  4Trr‘  MIS 
dose,  i.e., 
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K-Factor 


(10) 


2-D  4Trr^  variable  density  air  dose 
1-D  ^TTr^  MIS  dose 

Nonortho gonal  Expanding  Coordinate  Systen 

Since  atmospheric  density,  p(z),  decreases  exponentially 
with  increasing  altitude  (z),  radiation  mean  free  path  in- 
creases exponentially  with  altitude.  To  maintain  a mesh 
interval  which  is  constant  in  a mean  free  path  sense,  AIRDIF 
employs  a nonorthogonal  coordinate  system  which  expands  with 
altitude.  AIRDIF' s coordinates  are  defined  by  the  equations 


X 

1 

= r 

(11) 

X 

= 

(12) 

2 

X 

3 

= z 

(13) 

where  r,  $,  and  z are  as 

shown  in  Fig . 1 . 

The  transformation 

from  the  AIRDIF  coordinate  system  to  cartesian  coordinates 
is  given  by 

X = r cos  (1M-) 

y = r e^'^^  sin  (15) 

z = z (16) 

where  H in  Eq  (m-)  and  (15)  is  the  mesh  expansion  parameter. 
This  parameter  is  varied  with  source  altitude  to  control  the 
rate  of  expansion  of  the  coordinate  system  to  coincide  as 
closely  as  possible  to  the  rate  of  density  decrease  of  the 
1962  U.S.  Standard  Atmosphere. 

that,  if  the  mesh  expansion  parameter  is  made 
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Note 


infinite  in  magnitude,  Eq.  (1^),  (15)?  and  (16)  reduce  to 
normal  cylindrical  coordinates.  This  property  enables  AIHD 
to  also  operate  in  1-D  infinite  homogeneous  air.  In  this 
mode,  the  density  function  f(z)  is  set  equal  to  the  density 
ratio  p(z  )/p  where  z is  the  source  altitude,  and  ’the 

o U o 

density  gradient  is  set  to  zero.  These  operations  are  all 
performed  automatically  by  the  code  when  the  1-D  mode  of 
operation  is  selected.  A pictorial  representation  of  the 
AIRDIE  expanding  coordinate  system  is  shov;n  in  Fig.  1. 

Finite  Differencing 

Eq  (8)  is  represented  in  terms  of  continuous  functions 
and  is  finite  differenced  in  AIRDIF  for  nijmerical  solution. 
The  region  of  problem  definition  is  represented  by  a grid 
network  as  depicted  in  Fig.  2.  By  approximating  the  deriv- 
atives with  first  and  second  central  differences  (Ref  e-:lS4), 
Eq  (8)  is  replaced  by  a finite  difference  equation  at  each 
mesh  point.  These  are  nine  point  finite  difference  equa- 
tions in  2-D  variable  density  air  because  of  the  nonorthcg- 
onal  coordinates.  They  reduce  to  five  point  equations  in 
homogeneous  air. 

Matrix  Equation  Iterative  Solution 

The  resulting  finite  difference  equations  can  be  written 
in  matrix  form  as 

A FS  = (17) 

where  A = a block  tridiagonal  matrix,  which  is  a f'unction 
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of  the  mesh,  atmospheric  model,  and  energy  group, 
= the  source  vector  for  group  g and  includes  the 
scatter  from  both  collided  and  uncollided  energy 
groups, 

= the  unknown  scalar  collided  fluence  for  group  g. 

Eq  (17)  is  solved  group  by  group,  starting  with  the 
highest  energy  group,  using  a block  iterative  method  kno’ACi 
as  successive  line  overrelaxation  (SLOR)  (Ref  5-199)-  To 
start  the  iterative  process,  a guess  must  be  made  for  the 
unknown  fluence  in  group  g.  In  AIRDIF,  the  initial  guess  for 
the  uncollided  fluence  in  group  1 (the  highest  energy  group 
containing  source  particles),  is  set  equal  to  twice  the 
value  of  the  ■'uncollided  fluence  in  this  group.  For  all  other 
groups,  the  iterative  process  is  initiated  by  using  the  solu- 
tion to  the  previous  group  as  the  initial  guess. 

For  each  group  calculation,  iteration  continues  luntil 
the  fluence  value  converges  to  the  correct  answer.  Conver- 
gence is  checked  by  the  equation 


p (p )_  p ( P-  1 ) 


■ (p ) 


.01 


(IS) 


where  F is  the  collided  fluence  at  each  mesh  point  and  p is 
the  present  iteration  number.  If  the  inequality  is  satisfied 
for  five  consecutive  iterations,  the  iteration  is  terminated 
and  an  additional  check  for  proper  convergence  is  made.  This 
check  is  made  by  comparing  the  product  of  the  A matrix  and 
the  group  fluence  vector  F^  to  the  group  scarce  vector  at 
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n ^ ^ 


all  mesh  points.  If  the  difference  is  greater  than  C.10, 
iterations  are  resumed.  This  process  is  continued  until 
both  convergence  checks  are  satisfied. 


Storage  Method 

The  matrix  A in  Eq  (1?)  is  block  tridiagonal.  As  a re- 
sult, most  of  the  elements  have  a value  of  zero.  To  take 
advantage  of  the  sparseness  of  the  tridiagonal  macrix,  spe- 
cial storage  procedures  are  used  in  AIHDI?  to  conserve 
computer  memory.  Also,  special  multiplication  schemes  are 
used  to  eliminate  unnecessary  multiplication  by  zero  ele- 
ments. The  following  discussion  summarizes  the  structure 
of  the  A matrix  and  storage  of  its  elements  within  AIETIF. 

Referring  to  Fig.  3,  the  n'zmber  of  entries  on  each  row 
is  equal  to  the  n'amber  of  nonzero  fluer.ce  points  in  the 
mesh.  This  number  (NPTS)  is  the  product  of  the  niumber  of 
horizontal  mesh  points  on  each  row  CMITQH)  and  the  n’umber  of 
altitude  rows  (NRCW).  Although  each  row  of  matrix  A 
contains  NETS  points,  there  are  at  most  nine  nonzero  elements 
on  each  row. 

If  the  finite  difference  equations  are  written  from 
left  to  right  starting  with  the  bottom  row  and  each  successive 
row,  then  the  matrix  A can  be  -.vritten  in  the  block  tridiagonal 
form  sho'wn  in  Fig.  4-.  Note  there  are  3 nonzero  matrices  per 
row.  The  zero  matrices  indicated  in  this  figure  do  not  have 
to  be  stored.  The  subscripts  in  Fig.  correspond  to  a given 
altitude  row. 


• t 


Mote:  Each  entry  in  this  figure,  including  zeros,  is  a Matrix. 


Fig.  5 shows  the  storage  setup  for  an  internal  xatrix 


element.  !Iote  that  although  each  of  these  matrices  is  of 
dimension  rTHCR*NKOR,  there  are  at  most  three  nonzero  elements 
on  each  row.  Therefore,  they  can  ce  stored  as  a 3*?TH0R  ar- 


ray. The  total  storage 


!d  for 


the  block  tridiagonal 


matrix  A is  then  (TIROV/)  (3)  ^R*3)*  If  the  entire  A matrix 

had  been  stored,  it  would  have  been  necessary  to  store  (liHOR* 
!TR0W)^  words.  Thus  a substantial  savings  in  computer  memory 
is  realized  by  taking  advantage  of  the  sparseness  and  storing 
the  A matrix  in  block  tridiagonal  form. 
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Note:  "NflOR"  is  the  number  of  horizontal  mesh  points  minus  one. 


III.  AIRDI?  CODE 


The  AIRDIE  code  is  described  in  this  section.  The 
discussion  includes  the  structure  of  the  code,  conputer 
operational  requirements,  capabilities,  applications  and 


limitations.  The  code  is  written 


Structure 

Modularity . The  AIRDIF  code  consists  of  b.3  -:odules. 

The  modules  are  the  calling  program,  39  subroutine  subpro- 
grams including  an  abort  routine,  and  3 f’uncticn  subprograms. 


The  struct’ure  of  the  orogram  is  shown  in  the  module 


ov 


diagram  (Figs.  6 and  7) • A brief  description  of  each  module 
is  given  in  Appendix  A.  The  calling  program  (AIRDIF)  initi- 
ates execution  of  the  program  by  calling  the  subroutine  FFiEC 
and  stops  execution  of  the  code  upon  regaining  control.  Pro 


gram  execution  can  also  be  terminated  by  the  subroutine  AROR 
if  an  error  is  discovered  by  any  subprogram  during  execution 
The  subroutine  EXEC  controls  sequencing  of  the  remaining  mod- 
ules . Variables  are  passed  and  returned  as  formal  para- 


meters through  call 


statements  between  subprograms 


anc 


calling  subprogram. 

Hierarchy.  In  a structured  modular  program  such  as 
AIRDIF,  the  flow  of  control  between  modules  is  first  from 
left  to  right  by  file  and  then  from  top  to  bottom  by  rank 
within  any  given  file.  The  sequencing  of  control  for  the 
AIRDIF  code  is  illustrated  in  Figs.  6 and  7*  Beginning  with 
the  left  file,  the  module  AIRDIF  initiates  execution  by 
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calling  the  subprogram  EXEC.  Flow  is  again  from  left  to 
right  by  file.  The  third  file  contains  the  modules  DATAIN, 
CALC,  and  RESOUT.  Within  the  file  flow  goes  from  top  to 
bottom  by  rank.  Therefore,  DaTAIN  gains  control  next.  As 
flow  is  always  from  left  to  right  first  and  then  from  top 
to  bottom,  control  would  next  be  passed  DF:I,  etc.  The  exact 
sequence  of  modular  flow  control  by  standard  hierarchy  is 
given  by  the  underlined  numbers  below  each  of  the  modules  in 
Fig.  6. 

Versatility.  The  modularity  of  a structured  code  lends 
great  versatility  to  the  program.  Since  each  subprogram  is 
a separate  module,  any  module  can  be  replaced  or  modified 
without  affecting  the  operation  of  the  remaining  modules. 

For  example,  AIRDIF  in  its  present  form  employs  a special  ■ 
form  of  diffusion  theory  for  the  basic  calculations  of  the 


iterative  solution  for  a number  of  relatively  small  matrix 
equations  to  the  known,  exact  solution.  The  subroutine  DATAI.’I, 
which  reads  in  the  problem  data,  and  its  associated  subroutines 
were  checked  independently  for  their  ability  to  read,  process, 
and  output  data.  The  subroutine  'AESH  and  its  associated 
routines  were  checked  for  their  ability  to  collapse  the 
coordinate  system  to  normal  cylindrical  coordinates  for  one 
dim.ensicnal  hom.ogenecus  air  calculations.  The  performance 
of  checks  such  as  these  not  only  eliminates  problem.s  during 
code  developm.ent , but  also  enhances  the  credibility  of  the 
results  when  the  code  becomes  ooerational. 


’cm.Duter  Recuirem;ents 


The  ma.jor  determ.ining  factors  in  the  coriputati; 


01  ccm.- 


puter  costs  are  execution  tim.e  and  mem.ory  (storage)  require- 
m.ents.  AIRDI?  has  been  run  on  both  the  CDC  7500  and  the  CDC 
6630  ccm.puter  system.s  at  Kirtland  ARE,  :i.M.  and  the  CDC  6500 
at  Wright-Eatterscn  ARE,  Chio. 

Time  (CPU)  . Execution  tim.es  vary  depending  on  the  prob- 
lem.. The  m.ajcr  factors  causing  variation  in  the  tim.e  are  the 
num.ber  of  m.esh  poi.nts,  the  num.ber  of  energy  groups,  and  the 
rate  of  convergence  of  the  iterative  solution.  Typical  times 
for  execution  o.n  the  CDC  7600  are  on  the  order  of  60  seconds 
for  1500  m.esh  points  and  53  energy  groups.  Execution  on  the 
CDC  6600  typically  takes  four  to  five  tim.es  longer  than  the 
CDC  7603. 

Storage . The  mem.ory  requirer.ent  for  the  CDC  6600  is 


/ - ' .-h’  ’ 


r 


4. 

approxiaately  200  K OCTAL'.  This  provides  adequate  storage 
for  1500  mesh  points  and  61  particle  energy  groups.  A maxi- 
mum of  Vo  neutron  groups  and  21  gamma  groups  can  be  used. 

The  AIRDIF  program  requires  15  magnetic  disc  files 
for  local  storage  of  calculated  quantities  and  output.  The 


use  of  these 

files  allows 

multiple  use  of  many 

of  the 

large 

arrays  during 

various  step 

s of  operation  of  the 

calcul 

aticnal 

algorithms . 

This  greatly 

reduces  the  computer 

memory’ 

require- 

ments.  The  data  output  (storage  file  #21)  can  be  made  into 
a permanent  file  for  later  access  of  information  by  edit 
codes.  This  is  done  through  the  use  of  the  proper  control 
cards,  which  on  most  systems  consist  of  a permanent  fils  re- 
quest card  before  compilation  and  a catalog  card  after  execu- 


tion. 

In  addition, 

S- 

" 1 " mus  t 

be  specifi 

ed  on  the  third 

entry 

of  data  card 

(see  the 

data  card 

instructions  in  Sec 

IV). 

A list  of  the 

usage  of  each  of  the 

sented  in  Table  I. 

Associated  Edit  Programs 

: Several  data  editing  codes  v/ere  developed  to  access  and 

‘ use  the  data  generated  by  AIRDIF.  The  functions  of  these 

programs  could  have  been  included  as  modules  in  the  main 
program  but  were  made  separate  entities  to  reduce  computer 
time  and  storage  requirements  of  the  main  program.  This,  of 
course,  results  in  savings  of  real  turn  aro'und  time  for  the 
user  by  increasing  the  job  priority  on  time  sharing  systems. 

t The  GDC  7600  requires  220  K OCTAL  LCM  and  33  X OCTAL 
‘ SCM  for  5000  mesh  points  and  61  particle  energy  groups. 
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Disc  File  ^ 

5 

6 

10 
1 1 
12 

13 

14 

15 

1 6 

17  (temper ary) 
17  - Record  =1 
17  - Recori  -2 
13  - Record  =1 
13  - Record  #2 
19  - Record  #1 

19  - Record  #2 

20  - Record  ,#1 

20  - Record  #2 
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Table  I 

Disc  File  Usage 
Usaee 

Prograa  ir*?ut  file 

Prograci  cutout  file 

Cross  section  storage 

Direct  fluences  storage 

Slant  range  storage 

Horizontal  ranges  storage 

Mass  range  storage 

Collided  fluences  soorage 

Total  fluences  storage 

Group  doses  storage 

Neutron  doses  storage 

Ganna  doses  storage 

4irr^  neutron  doses  storage  (2-D) 

4ir’  ganna  doses  storage  (2-D) 

4-rrr^  neutron  doses  storage  (1-D) 
4ttt^  ganna  doses  storage  (1-D) 
Neutron  K-Factor  storage 
Ganna  X-Factor  storage 
Problen  paraneters,  total  fluences, 
and  neshing  infornation  for  one 
perxanent  data  file. 


Two  edit  codes,  MISFIT  a.cd  DOSCOMP,  were  written  to  be 
used  in  conj-'anction  with  AIRDIF.  Basic  theory,  descriptions, 
and  user  guides  for  these  codes  are  included  in  Appendices 
3 and  C.  Only  a brief  description  of  these  edit  codes  is 


orovided  here. 


The  MISFIT  edit  code  is  utilised  to  access  a 1-D  infi- 
nite homogeneous  air  data  file  produced  by  AIRDIF.  MISFIT 
first  uses  the  AIRDIF  data  to  compute  doses  and  then  least 
squares  fits  the  doses  to  Murphy's  radiation  transmission 
equation  (Ref  7:3c).  Murphy's  radiation  transmission  equa- 
tion is  given  by 

ln(T)  = A -A  X+A  X^+A  (X)’^‘+A  (X)'^’+A  (X)'^’^A  ln(X)  (19) 

1 2 3 4 5 6 7 

where  T is  the  ^+-r^  dose  at  mass  range  X,  and  Aj  thru  Ay  are 

the  coefficients.  These  coefficients  (MIS  coefficients)  are 
then  used  as  part  of  the  input  for  the  2-D  variable  density 
air  AIRDIF  run  (Ref  card  =10  and  11).  AIRDIF  utilizes  these 
coefficients  and  Murphy's  transmission  equation  to  compute 
the  1-D  MIS  dose  which  is  required  in  the  computation  of 
K-Factors,  Fig.  3 illustrates  the  sequence  for  the  use  of 
the  MISFIT  edit  program  with  AIRDIF  and  the  associated  input 
and  output  for  MISFIT. 

The  DOSCOMP  edit  program  also  accesses  a data  file  pro- 
duced by  an  AIRDIF  run  and  computes  doses  using  any  response 
f’unction  desired.  Additional  dose  calculations  and  K-factors 
(if  opted)  are  made  in  seconds  of  computer  time  rather  than 
minutes  that  would  be  recuired  for  another  AIRDIF  run.  The 


'I 
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time  savings  are  realized  because  DOSCOM?  does  not  have  to 
numerically  solve  the  diffusion  equation.  Utilization  of  the 
DOSCOMP  edit  code  with  AIRDI?  is  illustrated  for  dose  only 
computations  in  Pig.  9 and  for  dose  plus  K-Factor  computations 
in  Fig.  10. 

AIRDIF  Capabilities 

Some  AIRDIF  code  capabilities  have  already  been  mentioned 
briefly  in  the  introduction.  A more  detailed  description  of 
the  various  capabilities  of  AIRDIF  is  presented  here. 

Operational  Modes . AIRDIF  can  be  run  in  either  1 -D 
homogeneous  air  or  2-D  variable  density  air.  In  the  1-D  mode, 
the  nonorthogonal  expanding  coordinate  system  is  collapsed 
internally  to  normal  cylindrical  coordinates  by  setting  the 
mesh  expansion  parameter  (H)  in  Fqs  (14)  and  (15)  to  1 x 10^° 
km.  In  the  1-D  mode,  the  density  function  f(z)  is  set  equal 
to  a constant  which  is  equal  to  the  ratio  of  the  atmospheric 
densities  at  the  source  altitude  and  sea  level.  The  density 
gradient  g(z)  is  then  internally  set  to  zero. 

In  2-D  variable  density  air,  AIRDIF  can  be  used  to 
compute  the  radiation  environments  and  (if  opted)  the  K- 
Factors.  If  K-Factors  are  desired  a 1-D  air  run  must  be 
made  first,  and  the  1-D  data  must  be  stored  on  a permanent 
file.  MISFIT  is  then  used  to  get  the  MIS  coefficients  to 
Murphy's  transmission  equation.  As  mentioned  abovej  these 
coefficients  then  become  part  of  the  data  for  the  2-D  air 
run. 
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Problem  Parameters . AIRDIF  computes  the  free  field 


radiation  environments  for  neutrons,  prompt  gammas,  and 
secondary  gammas  for  atmospheric  nuclear  detonations.  Any 
source  spectra  can  be  used  for  the  computations.  The  only 
constraint  is  that  the  maximiim  number  of  groups  is  ^-0  for 
neutrons  and  21  for  gammas.  Doses  for  any  material  desired 
can  be  calculated  by  imputing  the  appropriate  response  func- 
tions in  the  proper  energy  group  structure.  Any  cross  sec- 
tion set  may  be  employed  as  long  as  they  are  input  in  the 
multigroup  format  described  in  Sec.  IV. 

Meshing  is  done  internally  by  AIHDIF.  The  only  para- 
meter that  must  be  specified  is  the  source  altitude.  In  1-D 
air  any  source  altitude  betveen  0 and  20  km  can  be  specified. 
In  2-D  variable  density  air,  the  source  height  has  a lower 
limit  of  about  300  meters  (AIRDIF  does  not  contain  an  air/ 
ground  interface)  and  an  upper  limit  of  20  km. 
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IV.  AIRDIF  User's  Guide 


This  section  presents  a user's  guide  for  the  AIRDIF 
code.  Subsections  are  included  on  input,  output,  and 
sample  problems. 

INPUT 

The  input  for  the  AIRDIF  code  consists  of  three  types 
of  cards.  These  are  the  control  cards,  the  program  deck,  and 
the  data  cards.  These  basic  units  are  ordered  as  shown  in 
Fig.  11. 

Control  Cards . The  control  cards  consist  of  cards 
giving  the  job  identification,  system  requirements,  control 
commands,  and  file  separators.  These  cards  will  not  be  dis- 
cussed because  they  are  unique  to  individual  computer  systems. 

Program  Deck.  The  program  deck  contains  all  the  program 
instructions.  As  pointed  out  earlier,  AIRDIF  is  written  in 
ANSI  Standard  Fortran.  The  program  deck  can  be  input  as 
a source  deck  or  as  a compiled  binary  deck.  Of  course,  if 
AIRDIF  is  stored  on  a permanent  file,  it  can  be  called  by 
control  cards  and  the  program  deck  would  not  be  required. 

Data  Cards . Data  cards  contain  the  information  neces- 
sary to  define  a particular  problem.  Problem  defining  para- 
meters, cross  sections,  SLOR  convergence  factors,  response 
functions , and  MIS  coefficients  must  be  input.  Detailed 
instrucpions  for  each  data  card  are  given  in  the  following 
paragraphs . 

1 . Card  #1 : This  contains  the  title  for  the  problem  in  a 
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1 2A6  format  (thus  the  title  is  entered  in  the  first  72 
columns  of  the  first  card) . 


2.  Card  #2:  This  card  contains  two  integer  entries  in  a 

213  format.  The  first  entry  is  the  number  of  neutron 
energy  groups,  and  the  second  entry  is  the  number  of 
gamma  energy  groups. 

3.  Card  #3:  This  card  has  tvo  entries  in  a E1 0.0,13  format. 

The  first  entry  is  the  source  height  (km)  in  El 0.0 
format.  The  second  entry  is  the  mode  of  operation  of 
the  code  which  may  be  either  1-D  homogeneous  air  or  2-0 
variable  density  air.  The  mode  is  specified  by  either 
entering  a "1”  or  a "2"  for  the  second  entry  of  this 
card  in  I3  format  (ie.,  in  column  13)- 

4.  Card  #4:  This  card  contains  two  entries  in  a 13,16 

format.  Entry  number  one  is  the  cross  section  print 
option  in  I3  format.  Entering  a blank  or  a zero  in 
column  3 results  in  the  following  action:  "Do  not 

print  cross  sections."  Entering  a "1"  results  in  the 
action:  "Bo  print  cross  sections."  Entry  number  two  in 
16  format  is  the  store  data  option.  This  entry  allows  fc: 
the  storage  of  data  on  a permanent  file  (Tape  21 ) that 
may  be  accessed  by  edit  codes  after  execution  is  com- 
pleted. Entering  a blank  or  a zero  in  column  9 results 
in  the  action:  "Do  not  create  a permanent  file." 

Entering  a "1"  in  column  9 results  in  the  action:  "Do 

create  a permanent  file."  If  a permanent  file  is  to  be 
created,  most  systems  require  a request  card  before 
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execution  and  a catalog  card  after  execution  in  the 
control  card  section. 

Card  #5'  Card  nuaber  five  is  a deck  of  cards  contain- 
ing overrelaxation  factors  required  by  the  3L0R  iter- 
ative solution  algorithm.  The  optimum  SLOR  factors 
for  each  energy  group  in  the  DLC-3I  structure  (Ref  3) 
were  experimentally  determined  (Ref  3)*  These  factors 
are  listed  in  Table  II  and  should  be  entered  in  6E10.0 
format  (i.e.  six  entries  per  data  card).  The  total 
number  of  entries  must  equal  the  total  number  of  par- 
ticle energy  groups  which  were  specified  in  card  #2 
(i.e.  member  of  neutron  groaps  + momber  of  gamma  groups). 
The  neutron  groups  are  listed  first  starting  with  the 
highest  energy  group.  The  first  gamma  group  follows  the 
last  neutron  group  on  the  same  card  'unless  that  card 
already  has  six  entries.  Each  line  of  Table  II  repre- 
sents the  entries  for  one  card  'using  the  DLC-31  struc- 
ture mentioned  above.  SLOR  factors  for  a different  en- 
ergy group  struct'ure  can  be  obtained  from  Table  II  by 
interpolation.  (See  Table  XVII,  Appendix  E,  for  TLT- 
energy  structure.) 

Card  #5:  Card  number  six  is  a deck  of  cards  containing 

the  source  spectra  in  a 6E13.0  format.  The  source  soectra 
for  the  neutrons  is  entered  first  (highest  enerry  grouc 
first),  followed  immediately  by  the  source  spectra  fer  the 
ga.m.ma  rays  (again  hig.hest  energy  group  first).  Aeain  entr" 
are  miade  at  the  rats  of  six  per  card  '/;ith  the  first  gamma 
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Table  II 

Group  SLOR  Factors 


1 .2 

1.2 

1 .2 

1 .2 

1.2 

1 .2 

1 .2 

1 .2 

1 .2 

1-3 

1 

1 .4 

1.2 

1 .0 

1 .0 

1 .0 

1 .2 

1 .2 

1.2 

1.2 

1.2 

1 .2 

1.2 

1 .2 

1.2 

1 .2 

1 .2 

1 .2 

1 .2 

1 .2 

1 .2 

1.2 

1 .2 

1.2 

1 .2 

1 .2 

1 .2 

1 .7 

1 .6 

1 .h 

1 .4 

1 .4 

1 .k 

1 .4 

1 .^■ 

1 .2 

1.2 

1 .2 

1.2 

1.2 

1.2 

1 .4 

1 .0 

1 .0 

1.0 

1 .0 

1 .0 

1.0 

roup 

follov/ing 

the  last 

neutron 

group  on  the 

sam' 

as  described  above.  The  total  nunber  of  entries  must  be 
equal  to  the  total  number  of  particle  groups  as  specifie 
on  card  #2.  Typical  veapon-like  energy  spectra  for  neu- 
trons and  prompt  gammas  in  the  DLC-31  structure  are 
included  in  Appendix  E. 

Card  #7:  On  this  card  the  title  for  the  dose  response 

function  type  is  entered  in  12A6  format.  A representa- 
tive set  of  tissue  and  silicon  response  functions  for 
neutrons  and  gamma  rays  is  given  in  Appendix  E. 

Card  #3:  This  is  a deck  of  cards  which  contains  the 

multigroup  neutron  response  functions.  The  numbers  are 
entered  in  6E10.0  format  (six  per  card).  Again,  the 
number  of  entries  must  equal  the  number  of  neutron  energ; 
groups  specified  on  card  #2.  The  responses  are  entered 
in  the  order  of  decreasing  energy,  i.e.  entry  is  the 
response  function  for  neutron  group  -1  which  is  the 
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highest  energy  group. 

9.  Card  #9'-  This  deck  of  cards  contains  the  gacma  re- 
sponse functions.  Entries  are  nade  in  the  sane  fornat 
(6E10.0)  and  order,  highest  energy  group  first,  as 
described  above  for  card  #8. 

10.  Card  #10:  The  seven  coefficients  to  the  fit  of  the 

neutron  dose  1-D  nass  integral  scaled  (MIS)  data  are 
entered  on  this  card  in  7S10.0  fornat.  The  coefficients 
as  discussed  in  Section  II  are  obtained  fron  the  edit 
program  MISFIT  which  performs  a least  squares  fit  to 
the  actual  data  from  a 1-D  homogenous  air  AIRDIF  run 

to  Murphy's  radiation  transmission  equation  (see 
Equation  19)* 

11.  Card  #11:  The  seven  coefficients  to  the  fit  of  the 

gamma  dose  1-D  MIS  data  are  entered  on  this  card  in 
7E10.0  format.  Coefficients  on  cards  10  and  11  are 
used  by  AIRDIF  only  during  2-D  air  calculations  where 
K-Factors  are  desired.  During  AIRDIF  r-ans  in  which 
K-Factors  are  not  desired  these  cards  must  still  be 
included  in  the  input  deck,  but  blank  cards  may  be  used. 

12.  Card  ?12 : This  card  contains  the  title  of  the  cross 

section  set  in  ITAf  f orr.at . 

13.  Card  #13:  This  card  is  a deck  of  cards  containing  the 

multigroup  cross  sections  for  both  neutrons  and  gammas 
which  are  read  in  at  sea  level  density  (Po=1.225F-3  g/ 
cmM . The  cross  sections  are  entered  in  a format  which 
is  a modified  version  of  FIDO  format  (Ref  7)  given  by 


FORMAT  (6(I2,A1 ,16,13), 5X, 13) 

This  format  calls  for  six  consecutive  cross  sections 
to  be  entered  on  each  card  and  entered  in  (12 , A1 , 16 , I3 ) 
format  (12  columns  each).  The  six  cross  section  entries 
on  each  card  are  followed  by  five  blanks  (5X)  in  columns 
73  through  77 • The  last  three  columns  of  each  card 
contain  the  card  number  (I3)  of  the  cross  section  deck. 
The  format  of  each  of  the  cross  section  entries  (12, A1, 
16,13)  consists  of  four  separate  entries.  The  first 
entry  (12)  designates  the  number  of  times  this  cross 
section  shall  be  repeated.  The  second  entry  (A1 ) is 
the  repeat/terminate  option  parameter.  If  this  entry 
is  a ”T",  cross  section  reading  is  terminated  immediate- 
ly. The  previous  cross  section  read  would  be  the  last 
in  the  set  stored  in  the  AIRDI?  program.  If  this 
repeat/terminate  option  parameter  is  an  "R”,  then  this 
cross  section  will  be  entered  the  number  of  times  des- 
ignated in  the  first  entry  of  the  cross  section.  If 
this  entry  is  a blank,  then  the  cross  section  will  be 
entered  in  program  storage  only  once,  the  number  in 
the  first  entry  is  ignored,  and  the  program  will  ad- 
vance to  read  the  next  cross  section.  The  third  number 
in  the  cross  section  entry  is  the  number  base  of  the 
cross  section  entered  in  I6  format.  The  fourth  element 
of  the  cross  section  entry  is  the  exponent  and  exponent 
sign  of  the  cross  section  entered  in  I3  format.  The 
cross  sections  are  read  group  by  group  starting  with 
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the  highest  energy  neutron  group,  proceeding  through 
the  neutron  groups,  and  finally  proceeding  on  through 


the  ganna  groups,  again  starting  with  the  highest 
energy.  The  order  in  which  the  cross  sections  must 
be  entered  for  each  energy  group  is  given  in  Table  III 


Table  III 

Order  of  Cross  Sections  for  FIDO  Format 


Position 

1 

2 

3 

4 

5 

6 
7 

etc . 


Cross  Section 
removal 
transport 
total 

scatter  (g-^g) 
scatter(g-l-^g) 
scatter  (g-2-t-g) 

scatter  (g-3-^g) 
etc . 


note:  ’’g’’  represents  the  energy  group  (#1  is  hi 

Samples  of  the  card  entries  for  each  of  the  13 
data  cards  listed  in  the  above  data  card  instruction 
are  presented  in  the  sample  problems  in  Appendix  D. 
contains  a summary  of  the  AIRDIF  input  data  cards. 


;). 


input 

list 

Table 


Outou~ 

The  AIRDIF  output  consists  of  a cross  section  listing  ( 
opted),  problem  definition  parameters  listing,  input  data 
listings,  mesh  information,  and  results  in  terms  of  particl 
doses  and  K-Factors. 

Cross  Section  Listing . A complete  listing  of  the  cros 
sections  used  for  calculations  can  be  printed  in  the  output 


Table  lY 


Suninary  of 

Input  Data  Cards  for  AIRDIF 

Card  # 

Format 

Descriotior. 

1 

12A6 

Title  card  for  program 

2 

213 

Nunber  of  groups  (neutrons/g 

3 

E10.0,I3 

Source  height,  dinension 

13,16 

Option  card 

5 

6E10.0 

SLOR  convergence  factors 

6 

6E10.0 

Source  specora  (neutrcns/gan 

7 

12A6 

Response  f’unctions  title 

8 

6E10.0 

Neutron  response  f’anctions 

9 

6E10.0 

Gaairia  response  f-unctions 

10 

7E10.0 

Neutron  MIS  coefficients 

1 1 

7E10.0 

Gamta  MIS  coefficients 

12 

12A6 

Cross  section  title  card 

13 

FIDO 

Cross  sections 
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if  the  proper  entry  is  made  on  the  option  card  (card  #4)  as 
discussed  in  the  instructions  for  data  cards.  The  DLC-31 
cross  sections  as  printed  in  this  fornat  are  given  in  Ap- 
pendix E.  Identification  of  the  cross  sections  for  any  col- 
umn (energy  group)  in  this  table  (Table  car.  be  r.ade 
by  referring  bo  Table  III. 

Problem  Definition  Parameters  Listine.  The  ncr.tabular 


information  unique  to  a particular  problem  is  contained  in 
this  section  of  the  data  listing.  First  is  the  problem  title 
exactly  as  entered  on  data  card  #1 • The  title  appears  in  the 
output  immediately  after  the  cross  section  listing.  The  re- 
maining problem  definition  data  appear  next  in  the  folloving 
order:  source  height,  number  of  neutron  groups  , rrumber  of  gam- 
ma groups,  the  grid  increment  given  in  terms  of  mass  range, 
the  maximum  horisontal  range  in  terms  of  mass  range,  the  max- 
imum vertical  range  from  the  so'orce  in  terms  of  mass  range, 
and  the  desired  dimension  of  the  problem  (1-D  cr  2-D).  Ac- 
tual examples  of  this  output  can  be  fo’und  in  the  output  for 
the  sample  problems  presented  in  Appendix  D. 

Incut  Data  Listings . The  next  portion  of  the  output  is 
a listing  of  all  the  tab’uLated  input  data.  The  items  includ- 
ed are  group  SLOP,  factors,  neutron  source  spectra,  gam.ma 
source  spectra,  neutron  response  f’jnctions,  gamma  response 
f’onctions,  neutron  MIS  coefficients,  gamma  MIS  coefficients, 
and  the  title  for  the  cross  sections  used.  Actual  examples 
of  this  output  can  be  fo'und  with  the  sample  problems  presented 
in  Appendix  D. 
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Mesh  Ir.f’or.'nati 


The 


nex 


portion  of  the  output  is  a 


susimary  of  the  nesh  as  calculated  by  the  AIHDIF  program. 

The  following  are  listed:  the  nesh  expansion  paraneter,  nesh 
dimensions  in  terms  of  number  of  horizontal  and  vertical 
points,  a s’nmmary  of  the  horizontal  mesh  at  sea  level,  a 
summary  of  the  vertical  mesh,  and  the  increments  (delta  Z) . 

Results . The  last  section  of  printed  output  is  the 
section  containing  the  results.  This  section  titled  ''Doses 
and  K-Factors"  lists  all  mesh  points  of  the  2-D  mesh  scheme, 
line  by  line,  beginning  with  the  bottom  line.  The  points  on 
eac’n  line  are  given  in  order  from  the  axis  of  symmetry  (Z- 
axis)  to  the  outer  boundary.  The  horizontal  point  n'umbers 


are 

listed 

'und 

a r» 

column  ”1'' 

and 

the 

vertical  1 

ine  numbers 

are 

listed 

uni 

■5  r 

the  col'umn 

!T  T 
•J 

in 

the  output . 

For  each 

poin 

,t , the 

out 

listing  gi 

ves 

the 

following 

information: 

altitude,  horizontal  range,  slant  range  from  the  so'urce,  mass 
range  from  the  so-orce,  doses,  doses,  and  K-?actors  (for 

both  neutron  and  gamma  particles). 

There  are  three  possible  types  of  data  which  can  be 
output  in  the  X-Factor  column.  V.'hen  the  AIRDIF  code  is  oper- 
ated in  the  1-D  air  mode,  K-Factors  are  not-  computed  and  all 
entries  under  the  K-Factor  column  are  appropriately  printed 
as  unity  (1.00E+C0).  As  pointed  out  earlier,  in  the  2-D  var- 
iable density  air  mode,  AIRDIF  can  be  used  to  compute  envi- 
rorznents  alone  or  environments  plus  K-Factors.  If  K-Factors 
are  required,  the  MIS  coefficients  must  be  included  in  the 
input.  If  the  code  is  being  used  to  produce  environments 
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only,  then  the  MIS  coefficients  should  be  input  as  zeros 
(blank  cards).  This  action  will  result  in  all  1-D  MIS  doses 
being  set  equal  to  unity,  and  hence,  the  K-?actors  will  assur.e 
a pseudo  value  equal  to  the  2-D  doses. 

Another  possible  AIRDI?  output  is  the  error  diagnostic 
statement.  If  the  program  detects  an  error  at  any  place, 
execution  is  terminated  through  the  subroutine  ABORT.  Rc 
further  normal  output  data  are  printed  after  the  point  of 
error  detection.  However,  the  output  does  include  a state- 
ment of  error  diagnosis  which  explains  the  error  and  gives 
the  subroutine  location  of  the  error. 

Sample  Problems 

To  further  illustrate  how  to  use  the  AIRDIP  code  and  the 
associated  edit  programs,  detailed  sample  problems  are  pre- 
sented in  Appendix  0.  Both  the  input  cards  and  some  repre- 
sentative portions  of  the  printed  output  are  listed  for  each 
of  the  problems.  Sample  problems  are  included  for  both  1-D 
homogeneous  air  and  2-D  variable  density  air  operation. 

Table  V summarizes  the  problems  included  in  Appendix  D. 
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Table  V 


Sample  Frcblems 


GVPLU ; 1 . 


2 . 

3. 

4. 


FIND: 


SOLUTION; 


2. 


7 


A fission  detonation  at  5 km  altitude  and  a 
therr.onuclear  detonation  at  15  kn  altitude. 

Oak  Ridge  unclassified  normalized  neuoron  sourc 
spectra . 

DLC-31  (37/21)  coupled  neutron/garrjT.a  cross  sec- 
tions . 

Silicon  response  functions. 

Free  field  doses,  4-:r’  doses,  and  K-Factcrs 
(1-D  and  2-D  air). 

Follow  the  user's  guide  to  enter  data  on  cards. 
Using  the  fission  source,  perform  a 1-D  AIF.DIF 
run  obtaining  1-D  results  a.nd  create  a oermanen 
file  of  the  results. 

Exercise  the  NISFIT  program,  to  obtain  MIS 
coefficients . 

Using  the  fission  source,  perform,  a 2-D  AI.RDIF 
run  obtaining  2-D  results  and  K-Factors. 

Using  the  therm.onuclear  source,  perform  a 2-D 
AIF.DIF  run  obtaining  environ.m.ents  only. 


1 


V.  Limitations 


Although  diffusion  theory  gives  extremely  good  results 
at  very  little  cost,  as  verified  by  Shulstad  (Ref  1),  there 
are  certain  irdierent  limitations  to  be  aware  of  when  using 
the  code.  This  section  briefly  discusses  some  of  these 
limitations . 


Boundary  Value  Tesradation 

Because  the  angular  dependence  of  the  fluence  in  dif- 
fusion theory  is  limited  to  linear  variations,  the  diffusion 
results  become  suspect  and  can  not  be  used  with  confidence 
within  a few  mean  free  paths  (20-30  g/cm^)  of  a boundary:  the 
conditions  that  are  applied  here  (zero  fluence  on  the  sides 
and  zero  partial  return  current  at  the  top  and  bottom)  are 
only  approximations  made  to  close  the  solution.  Therefore, 
the  results  on  the  first  two  rows,  last  tv/c  rows,  and  the 
last  four  points  on  any  row  should  not  be  used. 

.A-ltltuie  Limits 

In  1-T  air,  AIRTI?  may  be  r'un  at  any  altitude  between 
0 and  20  km.  In  2-T  air,  restrictions  are  imposed  on  the 


lowest  and  highest  source  alt 


^ -•  1 


" L w ^ 


0 ot 


contain  an  air/ground  interface. 


down 


■«'r. 


om  tne 


source  height,  AIRDI?  will  set  the  lover  altitude  boundary 
at  250  g/cm^  below  the  so-urce  altitude  but  never  lower  than 
50  g/cm^  above  gro’und  (Z=0).  The  user  should  restrict  the 
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source  altitude  during  2-D  runs  to  altitudes  greater  than 
300  meters.  The  upper  source  altitude  limit  is  20  km  (can 
be  equal  to  20  km) . This  limitation  arises  from  the  in- 
ability of  simple  diffusion  theory  to  treat  streaming  rad- 
iation which  characterizes  particle  transport  at  high  alti- 
tudes where  the  air  density  is  very  low. 
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Appendix  A 


Module  Descriotions 


This  appendix  contains  brief  descriptions  of  the  modules 
In  the  AIRDI?  program.  For  each  module,  the  purpose,  a 
brief  summary  of  the  calculations  performed  within,  and 
the  sequential  relationship  to  other  m.odules  Is  discussed. 

The  order  of  presentation  of  m.odule  descriptions  Is  governed 
by  the  hierarchy  depicted  In  the  m.odule  flow  diagram  In 
Figs.  6 and  7.  The  exception  Is  that  the  three  function 
subprogram.s  are  given  last. 

Program  AIRDIF 


AIRDIF  Is  the  m.aln  program,  which  Initiates  execution  of 
all  remaining  modules.  Through  the  progra.m.  statem.ent,  this 
m.odule  acquires  and  attaches  all  necessary  disc  files.  The 
calling  program,  passes  control  to  the  executive  routine, 
EXEC,  through  a call  statemient  and  term.lnates  the  progra.m. 
upon  return  of  control.  The  subprogram.  EXEC  is  the  only  ex- 
ternal subroutine  connected  directly  to  the  m.aln  program.. 
Subroutine  Subprogram.s 

EXEC.  This  is  the  executive  routine  for  the  AIRDIF 


progra.m. . EXEC  receives  its  control  from  the  main  program, 
and  transfers  control  as  necessary  for  data  acquisition, 
calculation,  or  proper  disposition  of  results.  This  is  ac- 
com.pllshed  by  calling  DATAIN,  CALC,  and  RESO'JT  respectively. 
EXEC  also  reserves  the  required  m.em.ory  for  the  entire  progra.m 
by  initializing  and  dimensioning  the  arrays  and  setting  their 
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naxirrari  si’es.  Data  input  cards  #1  and  ir2  are  also  read  ty 
this  routine. 

ABORT.  The  routine  ABORT  can  be  called  froT;  any  place 


in  the  progra™  "vhere  an  error  is  detected.  The  type  oT  error 
is  minted  uoon  detection  and  the  ciodule  •.•/hich  detected  the 


5 c u.^ f. c r 3 ns. "C 0 ( 


error  is  identified, 
by  .ABORT. 

DATAIR.  DATAIM  is  called  by  EXTl.  This  is  the  ccTtjr.and 
e fcr  initiating  the  reading  of  problem  definition  data, 


3Ut 


source  soectra,  response 


jns,  SLCR  convergence  f: 


and  cross  sections.  The  only  card  actually  read  by  this  sub- 


routine is  -12,  the  cross  seel 


-he  stor 


■ f cross 


sections  on  disc  fils  -13  also  takes  place  under  the  cottmand 


■-J, 


, and  DATC'JT. 


^his  routine  is  called  by  DATA-R.  Input  dat 


> a o V " n ! 


cards  =3  through  -'ll  are  read  by  this  routi.ne.  

ranges  in  terms  of  m.ass  units  (or),  the  horizontal  m.ass  range 
a,t  coaltitude,  and  the  vertical  range  uo  and  do'.vn  from,  the 
source  are  determ.ined  by  DR!!.  The  m.esh  interval  or  Tiesh 
ooint  increment  is  also  determ.ined  here.  The  only  external 
to  this  routine  is  ?CI!IT. 

FCI!iT.  This  is  a routine  called  by  DR!!  to  determine  the 
numher  of  .m.esh  points.  ?CI!IT  determ.ines  the  num.ber  of  points 
up  and  dovn'.  from,  the  source,  the  num.ber  of  horizontal  points, 
a.nd  the  total  num.ber  of  ooints.  A check  is  also  m.ade  to 
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assure  niaxi2i''air.  array  sizes  are  not  violated  by  the  nesning. 

CROSS . CROSS  is  called  by  DATAIN  to  store  the  cross 
sections  in  nevnory  in  the  order  shown  in  Table  II.  CROSS 
calls  the  external  FIDMOD. 

FIDMOD.  FIDMOD  is  called  by  CROSS  during  the  reading 
of  the  cross  sections  (Card  #13)  to  interpret  the  FIDO  fornat. 
FIDMOD  then  passes  the  cross  section  to  CROSS  for  storage. 
INCHES  is  the  only  external  called  by  FIDMOD. 

INCRES . The  routine  DICRES  is  called  by  FIDMOD  to  in- 
crease the  group  and  position  counters  to  ins-ore  storage  of 
cross  section  in  the  correct  storage  locations. 

DATOUT.  DATOUT  is  the  routine  called  by  DATAIN  to  print 
the  infornation  read  from  the  input  cards.  The  meshing  in- 
formation is  also  printed. 

CALC . CALC,  called  by  EXEC,  is  the  command  and  control 
routine  for  performing  calcolations . All  group  calculations 
are  performed  within  CALC.  For  a s-ommary  of  the  theory, 
see  Sec  II  of  this  report.  CALC  passes  control  to  other 
subroutines  for  many  of  the  specialized  computations.  Nine 
externals  are  called  by  CALC:  MESH,  GSTART,  VIRGIN,  SCAT, 


TOTAL,  DOSE,  R2D0SE,  MIS,  and  KFACT. 

MESH.  The  routine  MESH,  called  by  CALC,  is  the  execu- 
tive routine  for  determining  the  mesh.  Calculations  are  not 
performed  by  this  routine,  but  accc -nplished  by  calling  five 
externals:  VMESH  for  the  vertical  mesh,  SCALH  for  the  mesh 
expansion  parameter,  HMESK  for  the  horizontal  mesh,  RANGE 


for  the  distances  in  cm,  MASSRG  for  the  distances 


47 


•r-  wr' 


n terms 


of  mass  units.  Finally,  RZOUT  is  called  to  print  a mesh 
summary.  In  addition,  the  f'onction  subprogram  .A.IRDEN  is 
attached  by  VME3H,  SCALH,  and  I-IMESH.  The  subroutine  MESH 
also  insures  that  the  source  height  is  betv/een  0 and  20  km. 
Otherv/ise,  a diagnostic  m.essage  is  printed  and  control  passe; 
to  ABORT. 

GSTART.  The  routine  GSTAHT  is  called  by  CALC  to  deter- 
mine the  energy  group  to  start  calculations . This  routine 
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ergy  group.  The  computaticn  is  made  using  Eq . (1). 

SCAT . SC.AT  is  the  routine  called  by  C.ALC  to  determ.!: 
the  collided  fluence  at  each  mesh  point  of  the  entire  spa- 
tial domain  for  all  energy  groups  . T.ie  solution  is  dene 
group  by  group,  one  at  a time,  begiming  with  the  highest 


energy  group  containi.ng  particles 


relaxation.  SCA.  is  the  command  and  control  routine  for  the 
calculations  which  are  performed  by  calling  the  externals: 
DETDO,  MATRIX,  GUESS,  SCUHCE,  SLCR,  and  CHECK. 

DSTDO.  The  routine  LETDO  is  called  by  SCAT  to  compute 
the  group  diffusion  coefficients.  Eq.  (6)  is  used  to  perform 
the  computation. 

^-LATRIX . This  routine  is  called  by  SCAT  to  determine  the 
elements  of  the  block  tridiagonal  matrix.  These  are  the 
coefficients  of  the  finite  differenced  diffusion  equation. 
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The  elements  are  stored  in  the  order  discussed  in  Sec.  II  and 
illustrated  by  Figs.  3,  and  5*  liATRIX  calls  the  f'anction 
subprograms  AIRDEX  and  GINT. 

DENGRA.  This  routine  is  called  by  MATRIX  to  compute  the 
density  gradient  g(z)  at  altitude  (z)  for  each  mesh  row. 

GUESS . GUESS  is  the  routine  called  by  SCAT  to  make  the 
initial  guess  for  the  iterative  solution  of  each  energy  group. 
For  the  first  group  (highest  energy)  the  guess  is  set  equal 
to  twice  the  uncollided  fluence  for  the  corresponding  uncol- 
lided fluence  in  this  group.  For  other  groups,  the  answer 
for  the  previous  group  is  used  as  the  initial  guess. 

SOURCE.  This  is  the  routine  called  by  SCAT  to  compute 
the  soijnce  particles  for  each  energy  group  as  discussed  in 
Sec.  II.  SOURCE  calls  the  function  subprogram  AIRDEII. 

SLOR . SLCR  is  the  routine  called  by  SCAT  to  find  the 
solution  for  the  collided  fluence  at  each  mesh  point  for 
each  energy  group.  This  routine  uses  a block  iterative 
method  called  successive  line  overrelaxation  (SLOR)  to  solve 
the  matrix  equation  A*F  = S . SLOR  checks  for  convergence 
by  using  Eq  (13).  The  error  must  be  less  than  or  equal  to 
0.01  for  five  consecutive  iterations.  SLOR  calls  the  two 
routines  TRIDGE  and  TRIIFUL. 

TRIIiUL.  TRIIfJL  is  called  by  both  SLOR  and  3LTRIM.  This 
routine  is  a special  routine  for  multiplying  a tridiagonal 
matrix  by  a vector. 

TRIDGE . The  subroutine  TRIDGE  is  called  by  SLOR  and 
contains  a special  form  of  Gaussian  elimination  for  equations 
containing  a tridiagonal  matrix. 


CHECK.  CHECK  is  called  by  SCAT  to  check  for  proper  con- 
vergence of  the  entire  natrix  equation.  The  check  is  cade  by 
conparlng  the  tvo  sides  of  the  equation  A E = S for  the  last 
iteration.  If  the  difference  is  greater  than  0.10,  additional 
iterations  are  required  by  SLOP..  CHECK  calls  the  subroutine 
BLTHIM. 

BLTRIM.  BLTHIM  is  called  by  CHECK.  This  routine  mul- 
tiplies a vector  by  a block  tridiagonal  banded  matrix.  The 
routine  TRIMKJL  is  called  by  BLTRIM. 

TOTAL.  The  subroutine  TOTAL  is  called  by  CALC  to  find 
the  total  fluences.  This  is  accomplished  by  addition  of  the 
colXidsd  3.n.d  \incolXid6d  fX''sJi8nc0s  3."t  93.c''i  ir-9sri  oOir*'t  lor  63.cir. 
energy  group. 

DOSE.  DOSS  is  called  by  CA.LC  to  compute  doses.  The 
total  fluences  at  each  mesh  point  are  multiplied  by  the 
proper  response  Sanctions  to  give  the  dose  for  each  group  at 
that  point.  The  total  dose  summed  over  all  groups  is  com- 
puted in  CALC. 

R2D0SE . This  routine  is  called  by  CALC  to  compute  the 
Virr^  doses.  The  computation  is  accomplished  by  multiplying 
total  doses  at  each  mesh  point  by  '-+Tr^  where  r is  the  slant 
range . 

MXS.  MIS  is  a routine  called  by  CA.LC  to  compute  the 
1-D  MIS  4T7r^  doses  using  Eq  (19)  and  the  coefficients  that 
were  input  (cards  #10  and  #11).  The  f^xncticn  subprogram 
CU3RT  is  called  by  MIS. 

KF.ACT.  The  routine  KEACT  is  called  by  C.ALC  to  comoute 


K-Factors.  This  is  done  by  computing  the  ratio  of  the  2-D 
variable  density  air  ^:rr^  doses  divided  by  the  1-D  MIS  ^vr^ 
doses . 

RESOUT . RESOUT  is  called  by  EXEC  to  perform  disposi- 
tion of  the  results  of  the  computations.  RESOUT  calls 
DOSOUT  to  print  the  doses  and  K-?actors.  If  the  store 
option  on  input  card  #4  is  exercised,  STORE  is  called  to 
write  group  fluences  and  mesh  information  on  a permanent 
file  for  future  editing  operations. 


Function  Subprograms 

PINT.  GINT  is  a function  subprogram  called  by  ^-LA.SSRG. 
GINT  performs  Gaussian/Legendre  quadrat’ore  integration  on  th 
AIRDEN  function  to  compute  mass  integrals  for  variable  densi 
ty  air  (Ref  4:153-159)  (Ref  9:196). 

AIRDEN.  AIRDEN  is  a function  subprogram  called  by 
YMESH,  SCALE , rCIESH,  yJiSSF.C,  GINT,  MATRIX.  DENGRA,  and  SOURC 
AIRDEN  calculations  the  relative  air  density  o(s)/p(0)  as  a 
function  of  geometric  altitude.  The  calculations  are  based 
on  the  1962  U.S.  Standard  Atmosphere  (Ref  2)  and  are  valid 
over  the  altitude  range  from  5 km  below  sea  level  through  70 
km  above  sea  level,  inclusive. 

CU3RT . The  function  subprogram  CUBRT  is  used  by  MIS  in 
computation  of  the  cube  root  of  a given  quantity. 


Fig.  13.  Module  Flov  “iagran:  for  M.ISFIT. 


User's  Guide 


This  section  provides  instructions  for  preparation  of 
the  input  cards  for  MISFIT  and  presents  exacples  of  a typical 


input  data  deck  and  a listing  of  the  output 


The  basic 


' ^ ' 


r 


structure  of  the  deck  is  the  same  as  for  the  AIRDIF  code 
shown  in  Fig.  11.  There  are  three  input  data  cards.  Instruc- 
tions are  given  for  each  card  individually  below. 

Card  ^1  . This  is  the  title  card  for  the  response  f'onc- 
tion  type.  The  entry  is  alphameric  in  1 2A6  format  (same  as 
card  =7  in  AIRDIF  input). 

Card  -2 . This  is  a deck  of  cards  •vhich  contains  the 
neutron  response  functions.  The  rrombers  are  entered  in 
6E10.0  format  (same  as  card  in  AIRDIF  input). 

Card  . This  is  a deck  of  cards  containing  the  gamma 
response  functions.  Entries  are  made  in  6E1C.0  format.  This 
card  is  identical  to  card  #9  in  the  AIRDIF  input. 

Table  VI  summarizes  the  input  data  cards  for  kISFIT.  and 
Table  VII  depicts  specific  examples  of  punched  card  entries. 


Card  ^ 
1 
2 
3 


Summary  o:l 

Format 

12A6 

6S10.0 

6E10.0 


Table  VI 
It  Data  Ci 

Descriotion 

Response  f'unctions  title 
Neutron  response  f-mctions 
Gamma  response  f'unctions 


The  output  consists  of  a summary  of  the  problem  defi.ni- 
tion  and  input  data,  a listing  of  the  fit  values  for  the 
4-rr^  doses  at  various  distances  in  terms  of  mass  range  (or), 
and  listing  of  the  MIS  coefficienos  for  neutrons  and  gammas. 
The  listing  of  MIS  coefficients  also  gives  information  on 
the  m.aximium.  deviation  and  the  standard  deviation  of  the 
fit  and  the  nur.ber  of  points  rejected.  A sam.ple  of  the 
printed  output  is  given  in  Table  VIII. 
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Note;  See  list  on  page  64  for  notes  applicable  to  Table  VII. 

This  MISFIT  run  accessed  the  permanent  file  created  In  sample  problem  #1,  Appendix 


VERT  ROWS  = 52 
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Appendix  C 


DOSCOMP  Edit  Fro  eraa 


The  prinarv  piorpose  of  the  DOSCOMP  program  is  to  compute 
doses  with  a minimum  expendit'ore  of  computer  resources.  The 
AIRDIE  program  is  still  needed  to  make  the  basic  computation 
of  group  fluences  at  the  mesh  points.  However,  doses  for  only 
one  kind  of  dose  respo.nse  can  be  obtained  from  a single  AIRDIF 
execution.  Doses  for  the  same  source  but  of  a second  material 
response  would  require  another  complete  execution  of  AIRDIE 
which  requires  approximately  *+  to  ? minutes  on  the  CDC  6600. 
DOSCOMP  eliminates  this  needless  expenditure  of  computer 
time  by  accessing  a permanent  file  containing  data  from  the 
first  AIRDIE  execution.  The  reduced  time  and  decreased 
storage  requirements  of  DOSCOMP  greatly  decrease  the  turn- 
around time  for  computers  on  a time  sharing  system.  Herein 
lies  the  strength  of  the  DOSCOMP  edit  program  - it  allows  the 
computation  of  additional  doses  and  K-Eactors  with  a sub- 
stantial decrease  in  computer  storage  and  time  requirements, 
and  a large  reduction  of  real  time. 

Descriotion 

The  DOSCOMP  program,  like  AIRDIE,  is  a fully  structured 
modular  program  utilizing  standard  hierarchy.  The  program 
consists  of  13  modules,  many  of  which  were  extracted  directly 
from  AIRDIE.  A module  flow  diagram  for  the  DOSCOMP  program 
is  presented  in  Fig.  1U-. 


There  are  three  basic  sections  to  the  program.  These 
are  the  input,  calciolation,  and  output  sections.  VJhen  con- 
trol is  passed  to  the  input  section,  data  from  input  cards 
and  a permanent  file  are  read  into  computer  storage.  Next, 
control  is  passed  to  the  calculation  sections.  The  group 
fluences  are  then  used  to  compute  doses,  4Trr^  doses,  and  K- 
Factors  if  opted.  Finally,  the  results  are  printed.  Due 
to  the  modular  construction  of  the  program,  plot  subroutines 
could  easily  be  added  to  the  output  section. 


User's  Guide 

The  basic  structure  of  the  complete  card  deck  is  the 
same  as  given  in  Fig.  11  for  the  AIP.DIF  program.  The  perm- 
anent file  with  AIRDIF  data  must  be  accessed  using  the  proper 
control  cards.  Instructions  are  given  below  for  the  input 
data  cards. 

Card  --^1  . This  is  the  title  card  for  the  response  func- 
tion type.  The  entry  is  alphameric  in  12A6  format  and  is 
identical  to  card  #7  of  the  AIRDIF  input. 

Card  #2.  This  is  a deck  of  cards  which  contains  the 
neutron  response  functions.  The  nrambers  are  entered  in 
6E10.0  format.  This  is  identical  to  card  #3  in  the  required 
AIRDIF  input. 

Card  This  deck  of  cards  contains  the  gamma  response 

functions.  Entries  are  made  in  6S10.0  format  exactly  as 
accomplished  for  card  #9  of  the  AIRDIF  input. 

Card  The  neutron  coefficients  for  Murphy's  trans- 

mission equation  are  entered  on  this  card  in  7E10.0  format. 
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Surmiary 

Table  IX 

of  DCSCOX?  Input  Data  Cards 

Card  ^ 

Fornat 

Descriotion 

1 

12A6 

Response  function  title 

2 

6E10.0 

Neutron  response  f-Jinctions 

3 

6E10.0 

Ganna  response  functions 

if 

7E10.0 

Neutron  MIS  coefficients 

5 

7E10.0 

Ganna  MIS  coefficients 

6 

13 

Print  option  card 

DOS COM? 

Table 

P'onched  In 

X 

.put  Data  Cards 

Card  = 

Eornat 

Pm? 

.ched  cards 

1 

12A6 

MURRHY/5A 

RTINE  SILICON  RESPONSE 

2 

6E10.0 

sane 

as  card  ffS-Problen  , 

3 

6E10.0 

sane 

as  card  #9-?toblen  #1  , 

7E10.0 

sane 

as  card  #10-?roblen  #1 

5 

7E10.0 

s^m6 

as  card  #11-?roblen  #1 

6 

13 

bb1 
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Output  for  these  problems  are  sho’vn  in  Tables  XII,  XIV, 
and  XVI.  Although  all  output  sections  are  i nc  ed,  only 

limited  portions  of  lengthy  sections  are  included. 

Problem 

The  problem  statement  in  terms  of  given  quantities  and 
desired  results  is  given  in  Fig.  15*  Solution  instructions 
are  also  given  in  Fig.  15*  The  entries  for  the  input  data 
cards  to  this  problem  are  given  in  Table  XI.  A sample  of 
the  output  from  this  problem  is  given  irr  Table  XII.  The 
cross  section  listing  is  deleted  from  Table  XII,  but  is 
presen'^ed  in  Table  XX  in  Appendix  E. 

Problem  -2 

The  problem  statement  in  terms  of  given  quantities  and 
desired  results  is  given  in  Fig.  16.  Instructions  for  the 
problem  solution  are  also  given  in  Fig.  16.  The  input  data 
card  entries  for  this  problem  are  given  in  Table  XIII.  The 
output  from  this  problem  is  shcvn  in  Table  XIV. 

Problem 

The  problem  statement  is  given  in  Fig.  17  in  terms  of 
given  quantities  and  desired  results.  The  steps  for  the  so- 
lution are  also  given  in  Fig.  17  for  problem  ^3-  The  input 
data  card  entries  for  this  problem  are  given  in  Table  XT. 

The  output  is  shown  in  Table  XVI. 
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AIRDIF  - Frcblei  ^2 

GIVEN:  1.  A fission  detonation  at  5 kn  altitude. 

2.  Oak  Ridge  unclassified  normalized  neu- 
tron source  spectra  (Table  XVIII.  Ap- 
pendix E) . 

3.  DLC-31  (37/21)  coupled  neutron/gamma 
cross  sections  (Table  !Cv,  Appendix  S) . 

4.  Silicon  response  functions  (Table  XIX, 
Appendix  E) . 

5.  The  permanent  file  created  from  AIRDIE 
Problem  #1 . 

FIND:  1 . Neutron  and  secondary  gamma  environ- 

ments for  2-D  variable  density  air, 

2.  Produce  X-Factors. 

SOLUTION:  1 . Perform  a MISFIT  r'un  using  the  giver- 
response  f-onctions  and  the  permanent 
file  from  Problem  ^1.  (See  Table  VIII) 

2.  Perform  the  2-D  .IIRDIF  run  following 
the  user’s  guide  to  enter  data  on  the 
input  cards. 

3.  Desired  results  are  printed. 


ifig.  16.  -A.InDI?  Problem  -2  - Variabl 
Factors . 
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12A6  same  as  card  ^^2,  problem  #1 

FIDO  same  as  card  ^13,  pi’oblem 
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Tlie  "K"  for  tlie  exponent  lias  been  ilropiieil  fioin  tlie  "DOSES  AND  KFACTOKS"  section  of  this  table. 


Appendix  E 
Sanple  Input 

The  energy  group  structure  used  in  the  sarr.ple  problems 
and  listed  in  Table  XVII  is  the  DLC-31  (37/21)  .T.ultigroup 
structure  (Ref  3).  Source  spectra  corresponding  oo  these 
energy  groupings  are  given  in  Table  XVIII  for  neutrons 
(Ref  10)  and  prcnpt  ganna  particles  (Ref  3).  The  response 
functions  used  are  listed  in  Table  XIX  for  both  neutrons  and 
gamas  (Refs  11,  12,  and  13).  Table  XX  is  a listing  of  the 
cross  sections  used  in  the  sanple  problems.  These  are  the 
diffusion  theory  cross  sections  from,  the  Oak  Ridge,  DMA 
approved,  DLC-31  (37/21)  coupled  neutron/gam.mia  m.ultigrcup 
cross  section  set  (Ref  7) . 


Group  Neutron  Knort’iy  (Mcv)  Group  Neutron  Kner^.y  (Mev)  Group  Gamma  Erier^’iy  (Mev) 


Groups  jjB  Uiroufj;h  C>8  ai'o  the  prompt  (j;amma  enei);y  (^I'oups.  (From  Kef  8) 


Neutron  and  Prompt  Gamma  Source  Spectra 
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